global SSDIMed "/disk/agedisk4/medicare.work/miller-DUA50377/proj_ssdi"
* Settings
version 16
do "$SSDIMed/scripts/_auxiliary/_project_settings.do"

cap log close 
log using $SSDIMed/scripts/logs/04e_marginal_function_calculations.log, text replace 
cap mkdir "$SSDIMed/results/estimates/marginal_function"

****STEPS:
**Step 1: bring together true and bootstrap results for alphas and betas 
**a row in the data is a true or bootstrap sample
**Step 2 compute true or bootstrap slope and intercept for the lines in Figure 6
**Step 3: collapse retaining true values, standard deviations, and 5th and 95th percentiles 

****COLUMN 1: baseline spending, COLUMN 2: baseline mortality
**Step 1: bring together true and bootstrap results for alphas and betas 
**a row in the data is a true or bootstrap sample

foreach column in "1" "2"{
if "`column'"=="1"{
use "$SSDIMed/results/estimates/x-age52xUR/bootstrap_incidence_spec02"
rename (intercept UR age52 age52xUR) (alpha alphaUR alpha52 alpha52XUR)
merge 1:1 sample using "$SSDIMed/results/estimates/x-age52xUR/bootstrap_spending_01.dta", nogen 
rename (intercept UR age52 age52xUR) (beta betaUR beta52 beta52XUR)
  }
if "`column'"=="2"{
use "$SSDIMed/results/estimates/x-age52xUR/bootstrap_incidence_spec02"
rename (intercept UR age52 age52xUR) (alpha alphaUR alpha52 alpha52XUR)
merge 1:1 sample using "$SSDIMed/results/estimates/x-age52xUR/bootstrap_spending_01.dta", nogen 
rename (intercept UR age52 age52xUR) (beta betaUR beta52 beta52XUR)
}
**Step 2 compute true or bootstrap slope and intercept for the lines in Figures 6 & 7
***benefit parameters
gen double m=2*beta52/alpha52 
gen double n=beta-(alpha*beta52/alpha52)
gen double mUR=2*(beta52+beta52XUR)/(alpha52+alpha52XUR) 
gen double nUR=beta+betaUR- (alpha+alphaUR)*(beta52+beta52XUR)/(alpha52+alpha52XUR) 
gen double mURminusm=2*(beta52+beta52XUR)/(alpha52+alpha52XUR) - 2*beta52/alpha52
gen double nURminusn=betaUR- (alpha+alphaUR)*(beta52+beta52XUR)/(alpha52+alpha52XUR) + (alpha*beta52/alpha52)
***elements of the figure
gen double x51=alpha 
gen double x51UR=alpha+alphaUR 
gen double x52=alpha+alpha52 
gen double x52UR=alpha+alpha52+alphaUR+alpha52XUR
gen double B_x51=n+m*x51 
gen double B_x51UR=nUR+mUR*x51UR 
gen double B_x52=n+m*x52
gen double B_x52UR=nUR+mUR*x52UR
***difference at x51
gen double BU_x51=nUR+mUR*x51
gen double diff_x51=BU_x51-B_x51
***set up CI 
gen double B_200=n+m*200
gen double BUR_200=nUR+mUR*200 
gen double B_800=n+m*800
gen double BUR_800=nUR+mUR*800 
***cost parameters
foreach DeltaC in -1000 -5000{
local posDeltaC=-1*`DeltaC'
if "`DeltaC'"=="-1000" local span51=1 
if "`DeltaC'"=="-1000" local span52=20
if "`DeltaC'"=="-5000" local span51=.2 
if "`DeltaC'"=="-5000" local span52=4
gen double mC_51`posDeltaC'=(-`DeltaC' -m*alpha -n +mUR *(alpha+alphaUR) +nUR)/alphaUR
gen double nC_51`posDeltaC'=m*alpha +n- mC_51`posDeltaC'*alpha
gen double alphaCF_51_C`posDeltaC' =(`DeltaC' +nC_51`posDeltaC'-n)/(m-mC_51`posDeltaC')
gen double mC_52`posDeltaC'=(-`DeltaC'  +mUR*(alpha+alphaUR+alpha52+alpha52XUR) -m*(alpha+alpha52) +nUR - n)/(alphaUR+alpha52XUR)
gen double nC_52`posDeltaC'=m*(alpha+alpha52)+n - mC_52`posDeltaC'*(alpha+alpha52)
gen double alphaCF_52_C`posDeltaC' =(`DeltaC' +nC_52`posDeltaC'-n)/(m-mC_52`posDeltaC')
gen double C_x51_lo`posDeltaC'=  nC_51`posDeltaC'         +mC_51`posDeltaC'*(x51-`span51')
gen double C_x51_hi`posDeltaC'=  nC_51`posDeltaC'         +mC_51`posDeltaC'*(x51+`span51')
gen double C_x51UR_lo`posDeltaC'=nC_51`posDeltaC'+`DeltaC'+mC_51`posDeltaC'*(x51UR -`span51')
gen double C_x51UR_hi`posDeltaC'=nC_51`posDeltaC'+`DeltaC'+mC_51`posDeltaC'*(x51UR +`span51')
gen double C_x52_lo`posDeltaC'=  nC_52`posDeltaC'         +mC_52`posDeltaC'*(x52-`span52')
gen double C_x52_hi`posDeltaC'=  nC_52`posDeltaC'         +mC_52`posDeltaC'*(x52+`span52')
gen double C_x52UR_lo`posDeltaC'=nC_52`posDeltaC'+`DeltaC'+mC_52`posDeltaC'*(x52UR-`span52')
gen double C_x52UR_hi`posDeltaC'=nC_52`posDeltaC'+`DeltaC'+mC_52`posDeltaC'*(x52UR+`span52')
}

**Step 3: collapse retaining true values, standard deviations, and 5th and 95th percentiles 
 
cap mkdir "$SSDIMed/results/temp"
cap mkdir "$SSDIMed/results/tables"
***true 
preserve 
keep in 1 
xpose, clear varname 
save "$SSDIMed/results/temp/tablecolumn`column'_value.dta", replace 
restore 
***sd of bootstrap values 
preserve 
drop in 1 
gcollapse (sd) *, fast 
xpose, clear varname 
save "$SSDIMed/results/temp/tablecolumn`column'_sd.dta", replace 
restore 
***percentiles of bootstrap values 
preserve 
drop in 1 
gcollapse (p5) *, fast 
xpose, clear varname 
save "$SSDIMed/results/temp/tablecolumn`column'_lb.dta", replace 
restore 
preserve 
drop in 1 
gcollapse (p95) *, fast 
xpose, clear varname 
save "$SSDIMed/results/temp/tablecolumn`column'_ub.dta", replace
restore 
clear 
gen stat=""
foreach var in value sd lb ub{
append using "$SSDIMed/results/temp/tablecolumn`column'_`var'.dta"
replace stat="`var'" if stat==""
}
reshape wide v1, i(_varname) j(stat) string 
rename v1* * 
rename _varname var 
drop if var=="sample"
rename (lb sd ub value) (lb`column' sd`column' ub`column' value`column')
save "$SSDIMed/results/estimates/marginal_function/tablecolumn`column'.dta", replace 
}




****final table output 
use "$SSDIMed/results/estimates/marginal_function/tablecolumn1.dta", clear 
merge 1:1 var using "$SSDIMed/results/estimates/marginal_function/tablecolumn2.dta", nogen 
***drop figure elements 
gen dropme=inlist(var, "x51", "x52",  "", "B_x51UR", "B_x52", "B_x52UR")
replace dropme=1 if inlist(var, "alpha", "alpha52", "alpha52XUR", "alphaUR", "B_200", "BUR_200")
replace dropme=1 if inlist(var, "beta", "beta52", "beta52XUR", "betaUR", "B_800", "BUR_800")
***keep only deltaC=-5000, except for Panel C stuff 
local posDeltaC 5000
replace dropme=1 if regexm(var, "1000")
replace dropme=0 if regexm(var, "5000")
replace dropme=0 if regexm(var, "alphaCF")
replace dropme=1 if regexm(var, "hi") | regexm(var, "lo")
drop if dropme==1
drop ub* lb* dropme 
*reshape wide value sd, i(var) j(column) string 
gen counter=1 if var=="m"
replace counter=2 if var=="n"
replace counter=3 if var=="mUR"
replace counter=4 if var=="nUR"
replace counter=5 if var=="mURminusm"
replace counter=6 if var=="nURminusn"
replace counter=6.25 if var=="BU_x51"
replace counter=6.5 if var=="B_x51"
replace counter=6.75 if var=="diff_x51"
replace counter=7 if var=="mC_51`posDeltaC'"
replace counter=8 if var=="nC_51`posDeltaC'"
replace counter=9 if var=="mC_52`posDeltaC'"
replace counter=10 if var=="nC_52`posDeltaC'"
replace counter=11 if var=="x51UR"
replace counter=12 if var=="alphaCF_51_C1000"
replace counter=13 if var=="alphaCF_51_C5000"
replace counter=15 if var=="x52UR"
replace counter=16 if var=="alphaCF_52_C1000"
replace counter=17 if var=="alphaCF_52_C5000"

sort counter 
drop counter 
order var value1 sd1 value2 sd2 
export excel $SSDIMed/results/tables/marginal_function.xls, firstrow(variables) replace 



